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1  Statement  of  the  problem  studied  and  Significance 

This  project  concerns  development  and  application  of  finite  difference  methods  using  Cartesian 
grids  for  differential  equations  with  complicated  geometry  or/and  with  interfaces. 

A  model  problem  is  the  Poisson  equation 

V-(/3Vt i)  =  f{x,y),  (1) 

with  the  following  two  important  settings: 

•  The  solution  domain  is  arbitrary  with  Neumann  or  mixed  boundary  conditions,  see  Fig  1  (a). 
It  is  well  known  that  the  finite  difference  methods  generally  will  lose  accuracy  for  arbitrary 
boundaries. 

•  Interface  problems  defined  on  a  rectangular  domain,  see  Fig  1  (b).  The  coefficient  /3{x,y) 
may  be  discontinuous  across  the  interface.  Whenever  we  have  situations  involving  two 
different  materials,  or  the  same  material  but  in  different  states,  the  physical  properties  will 
be  different  across  the  interface,  we  have  to  deal  with  interface  problems.  For  example,  the 
interface  between  oil  and  water,  the  interface  between  a  fluid  and  a  gas  bubble  are  examples 
of  different  materials.  The  ice-melting,  crystal  growth,  and  other  solidification  problems  are 
examples  of  interface  problems  of  the  same  materials  but  at  different  states. 

The  problems  we  are  interested  in  come  from  different  applications  and  they  have  one,  several, 
or  all  of  the  following  properties: 
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•  There  is  one  or  more  complicated  geometries,  boundaries  or  interfaces. 

•  The  coefficients  of  the  differential  equations  may  be  discontinuous  across  interfaces. 


•  The  solution  and  its  derivatives  may  be  discontinuous  across  interfaces. 


•  There  may  be  discontinuous  or  singular  sources,  or  dipoles  along  interfaces,  for  example 


/  —  fc  + 


X)dX, 


where  T  is  a  curve  in  two  dimensions,  and  a  surface  in  three  dimensions,  6  is  the  Dirac-delta 
function  which  is  not  a  standard  mathematical  function  and  is  only  defined  in  terms  of  the 
distribution  theory. 


•  The  interfaces  may  be  fixed  or  moving  and  may  develop  topological  changes  such  as  merging 
or  breaking. 


1.1  Why  Cartesian  grids? 

It  is  true  that  the  finite  element  method  may  be  used  to  solve  the  problems  of  the  interest. 

However  we  want  to  develop  finite  difference  methods  for  the  following  reasons: 

•  There  is  almost  no  cost  in  the  grid  generation.  This  is  very  significant  for  moving  bound¬ 
ary/interface  problems. 

•  To  get  sharp  interfaces  and  solutions,  we  want  to  use  the  infinity  norm  instead  of  energy 
norm  used  in  the  finite  element  method.  A  method  converges  with  high  order  accuracy  in 
the  energy  norm  may  not  converge  at  all  in  the  infinity  norm. 

•  We  want  to  take  advantage  of  many  efficient  packages/solvers  which  are  written  for  Cartesian 
grids. 

•  Recently,  the  level  set  method  has  been  very  successful  for  a  number  of  moving  inter¬ 
face/boundary  problems,  especially  for  the  problems  with  topological  changes  and  three 
dimensional  problems.  The  level  set  function  works  best  with  Cartesian  grids. 

•  Finite  difference  schemes  are  generally  easier  to  learn  and  implement  than  the  finite  element 
method. 


2  Summary  of  the  most  important  results 

2.1  Algorithms  development 

2.1.1  Immersed  finite  element  methods  for  interface  problems 

We  have  developed  two  different  finite  element  spaces  for  interface  problems  and  developed  the 
corresponding  Galerkin  finite  element  methods  for  elliptic  interface  problems  using  Cartesian 
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grids.  Some  theoretical  results  about  these  finite  element  methods  are  also  obtained.  We  believe 
that  we  are  the  first  to  develop  these  new  finite  element  spaces  over  Cartesian  grids.  With  the 
new  finite  element  spaces,  we  are  able  to  derive  stable  and  accurate  numerical  methods  with 
sharp  interfaces,  see  [12,  4],  In  the  non-conforming  finite  element  space,  a  basis  function  is 
defined  almost  exact  the  same  as  the  standard  piecewise  linear  finite  element  space  except  that 
we  enforce  the  natural  jump  conditions  across  the  interface.  This  finite  element  space  have  the 
following  properties: 

•  It  has  the  same  support  as  the  standard  piecewise  linear  finite  element  space.  Therefore  it 
is  very  easy  to  programming. 

•  A  global  basis  function  may  not  be  continuous  across  the  edges  whether  the  interface  cuts 
through.  Therefore  it  is  a  non-conforming  finite  element  space. 

•  Existence  and  uniqueness  of  the  basis  function  has  been  established. 

•  Theoretically,  for  a  solution  of  an  elliptic  interface  problem,  we  can  construct  an  interpo¬ 
lation  function  that  approximates  the  solution  to  second  order,  and  the  first  derivatives  of 
the  solution  to  first  order  within  each  triangle. 

•  The  method  is  second  order  in  the  L2  norm  and  the  energy  norm.  In  terms  of  the  L°°  norm, 
a  linear  to  super-linear  convergence  is  expected. 

•  The  finite  element  method  can  be  directly  applied  to  three  dimensional  problems. 

In  the  conforming  finite  element  space  for  interface  problems,  we  have  shown  that  it  is  impos¬ 
sible  to  construct  a  conforming  finite  element  space  using  polynomials  which  satisfies  the  natural 
jump  conditions  while  it  has  the  same  support  as  the  standard  piecewise  linear  basis  function 
using  piecewise  polynomials.  As  a  breakthrough,  we  came  to  an  idea  to  extend  the  support  of 
the  basis  function  to  one  more  triangle  along  the  direction  of  the  interface.  In  this  way,  a  piece- 
wise  linear  and  conforming  finite  element  space  can  be  constructed  for  interface  problems.  We 
think  this  is  a  major  advance  in  this  area.  The  conforming  finite  element  space  has  the  following 
properties: 

•  Existence  and  uniqueness  of  the  basis  function  has  been  established. 

•  The  finite  element  space  is  a  conforming  finite  element  space.  Therefore  most  of  the  standard 
theories  for  conforming  finite  element  spaces  apply  directly.  High  order  elements  can  be 
constructed  accordingly. 

•  Theoretically,  for  a  solution  of  an  elliptic  interface  problem,  we  can  construct  an  interpo¬ 
lation  function  in  the  finite  element  space  that  approximates  the  solution  to  second  order, 
and  the  first  derivatives  of  the  solution  to  first  order  within  each  triangle 

•  The  method  is  second  order  in  the  L2  norm  and  in  the  energy  norm.  Numerical  examples 
show  that  it  is  also  second  order  accurate  in  the  the  L°°  norm. 
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2.1.2  Maximum  preserving  scheme  for  elliptic  interface  problems  and  applications 

Some  new  difference  methods  have  been  developed  for  elliptic  interface  problems  that  involve 
discontinuous  coefficients,  singular  source  terms,  non-smooth  or  even  discontinuous  coefficients 
across  an  arbitrary  interfaces.  The  new  methods  are  based  on  Cartesian  grids  and  satisfy  the  sign 
property  that  guarantees  the  discrete  maximum  principle  using  an  optimization  approach.  The 
resulting  linear  system  of  equations  obtained  from  the  new  methods  is  diagonally  dominant  and 
its  symmetric  part  is  negative  definite.  The  convergence  of  the  new  methods  have  been  proved 
by  constructing  the  comparison  functions. 

With  this  work,  we  have  a  unique  second  order  finite  difference  method  using  Cartesian  grids 
for  elliptic  interface  problems  that  allows  variable  discontinuous  coefficients.  This  work  has  solved 
the  stability  problem  for  the  original  immersed  interface  method  and  provided  theoretical  foun¬ 
dation  for  the  immersed  interface  method.  We  believe  the  new  method  is  the  first  second  order 
method  (in  the  maximum  norm)  using  Cartesian  grids  for  variable  discontinuous  coefficients  that 
can  guarantee  the  discrete  maximum  principal  and  second  order  convergence,  see  [11]. 

The  method  has  been  coupled  with  a  multigrid  solver.  We  have  written  a  general  subroutine 
for  this  methods  with  the  interface  been  represented  by  a  set  of  control  points. 

The  generalization  of  the  method  to  three  dimensional  problems  is  described  in  [2,  3]. 

The  generalization  of  the  method  to  diffusion  and  advection  equations  with  a  new  multigrid 
solver  is  described  in  [1], 

2.1.3  New  Formulations  and  algorithm  for  interface  problems  in  polar  coordinates 

The  new  methods  are  based  on  a  formulation  that  transforms  the  interface  problem  with  a  non¬ 
smooth  or  discontinuous  solution  to  a  problem  with  a  smooth  solution.  The  new  formulation 
leads  to  a  simple  second  order  finite  difference  scheme  for  the  partial  differential  equation  and 
a  new  interpolation  scheme  for  the  normal  derivative  of  the  solution.  In  conjunction  with  the 
fast  immersed  interface  method,  a  fast  solver  has  been  developed  for  the  interface  problems  with 
piecewise  constant  but  discontinuous  coefficient  using  the  new  formulation  in  polar  coordinate 
system.  The  work  is  described  in  [14]. 

2.2  Applications  of  new  algorithms 

2.2.1  Level-set  function  approach  to  an  inverse  interface  problem  in  shape  identifi¬ 
cation 

A  model  problem  in  electrical  impedance  tomography  for  the  identification  of  unknown  shapes 
from  data  in  a  narrow  strip  along  the  boundary  of  the  domain  is  investigated.  The  representation 
of  the  shape  of  the  boundary  and  its  evolution  during  an  iterative  reconstruction  process  is 
achieved  by  the  level  set  method.  The  shape  derivatives  of  this  problem  involve  the  normal 
derivative  of  the  potential  along  the  unknown  boundary.  Hence  an  accurate  resolution  of  its 
derivatives  along  the  unknown  interface  is  essential.  It  is  obtained  by  the  immersed  interface 
method.  The  work  for  2D  and  3D  problems  are  explained  in  [6]  and  [3]  respectively. 
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2.2.2  Autophobic  spreading  of  drops 

In  collaboration  with  J.  K.  Hunter  of  UC  Davis  and  H.  K.  Zhao  of  Stanford  University,  we 
formulate  a  model  for  the  spreading  on  a  surface  of  a  drop  that  deposits  an  autophobic  monolayer 
of  surfactant.  We  present  numerical  solutions  of  the  model  equations  using  an  immersed  interface 
method  and  a  level  set  method,  see  [5]. 


2.2.3  Immersed  Interface  Method  for  Navier-Stokes  equations  with  moving  inter¬ 
face. 

We  have  derived  the  jump  conditions  for  Navier-Stokes  equations  with  an  interface  in  2D  and 
3D,  see  [9].  The  immersed  interface  method  for  the  incompressible  Navier-Stokes  equations  with 
singular  forces  along  one  or  several  interfaces  in  the  solution  domain  has  been  developed  in  [8].  The 
new  method  is  based  on  a  second-order  projection  method  with  modifications  only  at  grid  points 
near  or  on  the  interface.  From  the  derivation  of  the  new  method,  we  expect  fully  second-order 
accuracy  for  the  velocity  and  nearly  second-order  accuracy  for  the  pressure  in  the  mavimiim  norm 
including  those  grid  points  near  or  on  the  interface.  This  has  been  confirmed  in  our  numerical 
experiments.  Furthermore,  the  computed  solutions  are  sharp  across  the  interface. 

Using  the  fast  immersed  interface  method  and  the  vorticity  stream-function  formulation,  we 
have  developed  a  fast  method  for  Navier-Stokes  equations  defined  on  general  domains  in  [13]. 

2.3  Software  development 

We  have  written  some  public  subroutines  of  fast  solvers  for  Helmholtz  and  Poisson  equations  on 
irregular  domains  either  exterior  or  interior  with  various  boundary  conditions  have  been  developed 
and  made  public  through  anonymous  ftp. 

These  subroutines  use  the  fast  immersed  interface  method  by  introducing  a  unknown  jump 
either  in  the  solution  or  the  normal  derivative.  Using  the  fast  solver  from  Fishpack,  we  implicitly 
solve  a  linear  system  of  equations  using  GMRES  iteration.  These  methods  are  called  fast  solvers 
because  of  the  number  of  iterations  is  independent  of  mesh  size.  While  there  are  similar  meth¬ 
ods  for  problems  with  Dirichlet  boundary  conditions,  our  approach  for  problems  with  Neumann 
boundary  conditions  using  some  preconditioning  techniques  is  not  only  unique,  but  also  necessary 
to  guarantee  the  number  of  iterations  to  be  independent  of  mesh  size.  These  packages  include 


•  Fast  solves  for  Helmholtz/Poisson  equations  defined  in  an  arbitrary  interior  domain  with 
Dirichlet  boundary  condition. 

•  Fast  solves  for  Helmholtz/Poisson  equations  defined  in  an  arbitrary  exterior  domain  with 
Dirichlet  boundary  condition. 

•  Fast  solves  for  Poisson  equations  defined  in  an  arbitrary  exterior  domain  with  Dirich- 
let/Neumann  boundary  condition. 
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Our  method  also  works  for  Helmholtz  equations  defined  in  an  arbitrary  exterior  domain,  but  the 
number  of  iteration  depends  oh  mesh  size.  Whether  we  can  derive  a  fast  solver  for  this  type  of 
problems  is  an  open  question. 

Other  applications  and  contributions  can  be  found  in  [10,  15,  16,  7]  etc.. 
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Yaw  Kyei,  Guo  Chen,  Yan  Gong  (Ph.D  graduate  students  at  NCSU,  three  of  them  will  graduate 
in  the  summer,  2003). 


References 

[1]  L.  M.  Adams  and  Z.  Li.  The  immersed  interface/multigrid  methods  for  interface  problems. 
SIAM  J.  Sci.  Comput.,  24:463-479,  2002. 

[2]  S.  Deng.  Immersed  interface  method  for  three  dimensional  interface  problems  and  applica¬ 
tions,  Ph.D  thesis.  North  Carolina  State  University,  2000. 

[3]  S.  Deng,  K.  Ito,  and  Z.  Li.  Three  dimensional  elliptic  solvers  for  interface  problems  and 
applications.  J.  Comput.  Phys.,  184:215-243,  2003. 

[4]  R.  Ewing,  Z.  Li,  T.  Lin,  and  Y.  Lin.  The  immersed  finite  volume  element  method  for  the 
elliptic  interface  problems.  Mathematics  and  Computers  in  Simulation,  50:63-76,  1999, 

[5]  J.  Hunter,  Z,  Li,  and  H,  Zhao.  Autophobic  spreading  of  drops,.  J,  Comput.  Phys.,  183:335- 
366,  2002. 

[6]  K.  Ito,  K.  Kunisch,  and  Z.  Li.  Level-set  function  approach  to  an  inverse  interface  problem. 
Inverse  Problems,  17:1225-1242,  2001. 

[7]  K.  Ito  and  Z.  Li.  The  immersed  interface  method  for  nonlinear  interface  problems  with  an 
application,  submitted,  2002. 

[8]  M-C.  Lai  and  Z.  Li.  The  immersed  interface  method  for  the  Navier-Stokes  equations  with 
singular  forces.  J.  Comput.  Phys.,  171:822-842,  2001. 

[9]  M-C.  Lai  and  Z.  Li.  A  remark  on  jump  conditions  for  the  three-dimensional  Navier-Stokes 
equations  involving  an  immersed  moving  membrane.  Applied  Math.  Letters,  14:149-154, 
2001. 

[10]  Z.  Li  and  W.  Cai.  A  level  set-boundary  element  method  for  simulation  of  dynamic  powder 
consolidation  of  metals.  In  Numerical  Analysis  and  Its  Application,  edited  by  L.  Vulkov,  J. 
Wasniewski,  and  P.  Yalamov,  pages  527-534.  Springer- Verlag,  2001. 


7 


ff 


[11]  Z.  Li  and  K.  Ito.  Maximum  principle  preserving  schemes  for  interface  problems  with  discon¬ 
tinuous  coefficients.  SIAM  J.  Sci.  Comput,  23:1225-1242,  2001. 

[12]  Z.  Li,  T.  Lin,  and  X.  Wu,  New  Cartesian  grid  methods  for  interface  problem  using  finite 
element  formulation.  Numer.  Math.,  in  press,  also  NCSU  CRSC-TR99-5,  1999. 

[13]  Z.  Li  and  C.  Wang.  A  fast  finite  difference  method  for  solving  Navier-Stokes  equations  on 
irregular  domains.  J.  of  Commu.  in  Math.  Sci.,  1:180-196,  2003. 

[14]  Z.  Li,.  W-C.  Wang,  I-L.  Chern,  and  M-C.  Lai.  New  formulations  for  interface  problems  in 
polar  coordinates.  SIAM  J.  Sci.  Comput.,  in  press,  2002. 

[15]  S.  R.  Lubkin  and  Z.  Li.  Force  and  deformation  on  branching  rudiments:  Cleaving  between 
hypotheses.  Biomechanics  and  Modeling  in  mechanobiology ,  1:5-16,  2002. 

[16]  H.  Ly,  M.  Jolly,  F.  Reitich,  K.  Ito,  and  H.  T.  Banks.  Piecewise  linear  models  for  field- 
responsive  fluids.  IEEE  Transaction  on  Magnetics,  37:558-560,  2001. 


8 


